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Abstract 

Background: By reshuffling genomes, structural genomic reorganizations provide genetic variation on which 
natural selection can work. Understanding the mechanisms underlying this process has been a long-standing 
question in evolutionary biology. In this context, our purpose in this study is to characterize the genomic regions 
involved in structural rearrangements between human and macaque genomes and determine their influence on 
meiotic recombination as a way to explore the adaptive role of genome shuffling in mammalian evolution. 

Results: We first constructed a highly refined map of the structural rearrangements and evolutionary breakpoint regions 
in the human and rhesus macaque genomes based on orthologous genes and whole-genome sequence alignments. 
Using two different algorithms, we refined the genomic position of known rearrangements previously reported by 
cytogenetic approaches and described new putative micro-rearrangements (inversions and indels) in both genomes. A 
detailed analysis of the rhesus macaque genome showed that evolutionary breakpoints are in gene-rich regions, being 
enriched in GO terms related to immune system. We also identified defense-response genes within a chromosome 
inversion fixed in the macaque lineage, underlying the relevance of structural genomic changes in evolutionary and/or 
adaptation processes. IVloreover, by combining in silico and experimental approaches, we studied the recombination 
pattern of specific chromosomes that have suffered rearrangements between human and macaque lineages. 

Conclusions: Our data suggest that adaptive alleles - in this case, genes involved in the immune response - might 
have been favored by genome rearrangements in the macaque lineage. 

Keywords: Genome shuffling. Inversions, IVlacaque, Recombination, Adaptation, Meiosis, Tandem repeats. Evolutionary 
breakpoints 



Background 

Large-scale genomic changes, such as inversions, translo- 
cations, fusions and fissions, contribute to the reshuffling 
of the genomic architecture of organisms, providing new 
sources of variation on which natural selection can work. 
In recent years, there has been an increasing number of 
studies focusing on the role of chromosomal reorganiza- 
tions in adaptation and speciation processes [1-4], and 
more specifically on the influence of genome shuffling in 
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recombination ([5] and references therein). In this frame- 
work, the "suppressed recombination" model has provided 
compelling evidence and a theoretical framework to ex- 
plain how chromosome rearrangements are involved in 
speciation [6,7]. Under this model, reorganizations such as 
inversions would have a minimal influence on fitness 
when present in the heterokaryotype, but rather would 
suppress recombination between genomic regions in- 
volved in reorganization, leading to the reduction of gene 
flow between diverging populations. In this context, 
chromosomal rearrangements would act as genetic bar- 
riers, interfering in the fixation of favorable alleles and 
allowing for the accumulation of genetic incompatibilities 
[8]. As a way to test this hypothesis, subsequent studies 
have analyzed sequence divergence (patterns of nucleotide 
differentiation) between species as an indirect estimation 
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of recombination [9,10]. High rates of sequence divergence 
detected in genes located at, or near, chromosomal rear- 
rangements have often been interpreted as indirect evi- 
dence of chromosomal speciation through suppressed 
recombination [9-16]. However, few empirical data have fo- 
cused on the relationship between evolutionary breakpoint 
regions (EBRs) and recombination rates. Initial studies in 
Drosophila described a strong reduction of recombination 
around inversion breakpoints and within the reorganization 
itself [17]. Analogous studies in mammals are scarce, and 
the role of evolutionary regions in recombination has just 
started to be elucidated [5,14]. 

Whole-genome comparisons of distantly related mam- 
malian species have provided the basis for establishing 
models that can explain genome dynamics [18-22]. In this 
sense, different approaches have been developed in recent 
years to define homologous synteny blocks (HSBs; i.e., 
regions where gene order has been conserved among 
species) and EBRs (regions where the synteny has been 
disrupted by chromosomal reorganizations) among 
mammalian genomes [19,23-26]. Such reconstructions 
have revealed that genomic regions implicated in struc- 
tural changes which occurred during the evolution of 
species are not distributed randomly through the gen- 
ome, but instead they are clustered in regions that are 
more prone to break and reorganize [19,23-26]. The 
fact that some chromosomal regions have been reused 
during mammalian chromosomal evolution questions 
(i) whether these regions are physically labile due to 
their DNA sequence and/or structural chromatin con- 
formation, and (ii) whether they represent regions where 
selection against breakpoints is minimal [26]. Regarding 
the first assumption, previous studies on mammalian 
genomes have provided compelling evidence that EBRs 
can be linked to the presence of repetitive elements, 
such as transposable elements, segmental duplications 
and/or tandem repeats [19,25-31]. However, given the di- 
versity of repetitive elements in EBRs, it is likely that se- 
quence composition is not alone influencing genome 
instability, clamoring for the involvement of additional 
factors such as the state of the chromatin (i.e., open 
chromatin may drive chromosomal reorganizations [32]) 
or selective constraints. In this latter case, comparative 
genomic studies have shown that mammalian EBRs tend 
to localize in gene-dense regions [22,28,32]. But there is a 
long-standing debate on the mechanisms behind this well- 
known phenomenon. Several lines of evidences indicate 
that EBRs are precisely located between genes (i.e., 
intergenic regions, see [32]) not necessarily affecting 
gene structure/function, while others have reported 
possible gene expression changes due to genome 
reshuffling (see [33]). 

Given this context, the general picture of the genomic 
features and DNA organization of genomic regions 



affected by structural reorganizations is still incomplete, 
as is that of how genomic changes are transmitted to the 
offspring during the formation of germ cells and contrib- 
ute to speciation. If genomic shuffling does affect evolu- 
tionary processes through the mechanical shearing at 
evolutionary breakpoints, how does it impact on meiotic 
recombination? In this sense, the analysis of the most re- 
cent human and chimpanzee recombination maps has re- 
vealed that rearranged chromosomes presented lower 
recombination rates than chromosomes that did not suffer 
any reorganization since the human-chimpanzee common 
ancestor [5]. Elucidating upon whether this pattern also 
holds for other mammalian species would have a relevant 
impact on our understanding of the role of genome shuf- 
fling in speciation. Here, we have analyzed the effect of 
genomic structural changes on genetic recombination in 
the rhesus macaque to understand the mechanisms 
underlying chromosomal evolution in mammals and de- 
termine, in the long-term, the influence of chromosomal 
reorganizations on meiotic recombination. To this end, 
we have firstly characterized the genomic regions involved 
in chromosomal rearrangements between human and ma- 
caque genomes. The rhesus monkey {Macaca mulatta, 
Tribe Papionini, Catarrhini) is a primate species widely 
used in both biomedicine and evolutionary studies 
[34-37]. AU members of the Tribe Papionini {Macaca, 
Papio, Mandrillus and Cercocebus) are characterized by 
highly stable karyotypes that have been regarded to retain 
the ancestral Catarrhini karyotype [34,35]. Due to such 
characteristics, the macaque has often been used as a reli- 
able primate out-group candidate for evolutionary studies 
when studying great apes. But despite its importance, litde 
effort has been made in characterizing the genomic land- 
scape of HSBs and EBRs in this species since the initial re- 
lease of the rhesus macaque genome [32,36,37]. Here we 
provide a detailed genomic map of the structural rear- 
rangements between human and macaque. We have re- 
fined the genomic position of known rearrangements 
previously reported by cytogenetic approaches and de- 
scribed new putative micro-rearrangements (inversions 
and indels) between human and macaque genomes. More- 
over, we have analyzed the repetitive DNA content and 
gene density in relation to chromosomal reorganizations, 
as well as the effect of inversions in meiotic recombin- 
ation, detecting immune-related genes in evolutionary 
breakpoint regions in the macaque genome. 

Results and discussion 

Homologous synteny blocks (HSBs) and evolutionary 
breakpoint regions (EBRs) In human and rhesus macaque 
genomes 

In order to analyze the chromosomal reorganizations 
(fusions/fissions, translocations and inversions) between 
human and rhesus macaque, two different algorithms 
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were applied: SyntenyTracker and Cassis (See Material 
and Methods for further information). Both approaches 
can detect HSBs and EBRs in both species based on the 
positions of orthologous genes. By analyzing a total of 
16,133 orthologous genes between Homo sapiens (HSA) 
and Macaca mulatta (MMU), SyntenyTracker detected 59 
EBRs in both genomes, with a median length of 259 Kbp 
in the human genome (ranging from 11.6 Kbp to 4.6 
Mbp) and 163.6 kbp in the macaque (ranging from 16.9 
Kbp to 6.1 Mbp) (Table 1 and Figure 1). Cassis, on the 
other hand, detected 109 EBRs in the human genome, and 
111 EBRs in the macaque (Table 1), with a median length 
of 51.3 Kbp for human EBRs (ranging from 3 bp to 4.5 
Mbp), and 26 Kbp for the rhesus macaque (ranging from 
7 bp to 512.5 Kbp) (Table 1). Merging the results obtained 
with both approaches resulted in 121 EBRs in the human 
genome and 125 EBRs in the macaque. After applying a 
conservative filtering process (See material and methods 
for further information), we obtained a total of 74 EBRs 
and 101 HSBs in the human genome with a median EBR 
length of 51.6 Kbp (ranging from 3 bp to 508.7 Kbp). 
Likewise, we obtained 77 EBRs and 94 HSBs in the ma- 
caque genome (Table 1 and Additional file 1: Table SI, 
Additional file 2: Table S2 and Additional file 3: Figure 
SI), with a median EBR length of 26 Kbp (ranging from 
7 bp to 420.3 Kbp). 

Our own study represents a departure from those con- 
ducted previously [27,28,30] in that it relies on a detailed 
comparison between human and rhesus macaque genomes 
based on orthologous genes and whole-genome sequence 
alignments. Previous cytogenetic studies delineated the pri- 
mate ancestral karyotype, defining conserved syntenies 
among species and the direction of chromosomal rear- 
rangements in a phylogenetic context [24,34,38,39] . Species 
from the Tribe Papionini - including Macaca, Papio, Man- 
drillus and Cercocebus - are characterized by sharing the 
same karyotype and large-scale chromosomal reorganiza- 
tions since their divergence from a common primate an- 
cestor [34,35]. When comparing these species with the 
human karyotype, previous cytogenetic studies described 
the presence of 20 intra- and inter-chromosomal reorgani- 
zations [34,35,37,40]. These rearrangements include 12 
pericentric inversions affecting eleven chromosomes, four 
paracentric inversions involving four chromosomes and 
four fusions/fissions [34,35,37,40]. Overall, our in silico 



approach confirmed the presence of the above-mentioned 
macro-reorganizations and, thus, refined the breakpoints 
involved in both genomes (Additional file 1: Table SI). 
Moreover, we identified 39 and 41 previously undetected 
EBRs in the human and rhesus genome, respectively, 
affecting 13 different chromosomes (Additional file 1: 
Table SI). This resulted in 21 previously undetected reor- 
ganizations in the human genome and 23 in the macaque. 
Although previous cytogenetic studies have reported that 
sbc chromosomes (MMU6, MMU8, MMUll, MMU17 and 
MMU19) have been maintained coUinear between both 
species, our results suggest that only chromosomes MMU6, 
MMU12 and MMU17 have maintained a complete con- 
served synteny. Six new indels (insertions or deletions) were 
also identified, ranging from 3.3 Kbp to 2,784.8 Kbp in five 
different chromosomes (Additional file 1: Table SI). 



Tandem repeats do not accumulate in EBRs in the 
macaque genome; instead, they are correlated with the 
evolutionary history of chromosomes 

Previous comparative genomic studies have revealed that 
mammalian breakpoint regions are especially rich in re- 
petitive elements, such as segmental duplications [41-43], 
repetitive sequences [25,42], transposable elements and 
long regulatory regions [26,29,44-46]. We tested whether 
this pattern applies to the rhesus macaque, a species that 
has retained the ancestral Catarrhini karyotype [34,35]. To 
this end, the genome distribution of tandem repeats (TR) 
in this species was analyzed. A total of 701,128 loci, repre- 
senting 60.9 Mbp of the whole genome was detected. In 
order to study the genome-wide distribution of TR, the 
number of base pairs involved in TR by screening non- 
overlapping windows of 100 Kbp along the genome was 
analyzed. When considering the number of base-pairs 
(per 100 Kbp) of each window covered by TR, we 
observed a significant increase of TR in telomeric and 
centromeric regions when compared to HSBs and EBRs 
(Kruskal-Wallis test, p-value < 0.0001), mirroring previous 
observations obtained in the human and great apes 
[25,26]. Subsequently, we tested whether there was a 
correspondence between TR and EBRs in all macaque 
chromosomes. Our analysis indicated that EBRs are not 
significantly enriched by tandem repeats when compared 
to HSBs (Mann- Whitney U test, p-value > 0.05). 



Table 1 Summary of the evolutionary breakpoint regions (EBRs) detected in the human (HSA) and macaque (IVIIVIU) genomes 



SyntenyTracker 



Cassis 



Merged 



HSA 



MMU 



HSA 



MMU 



HSA 



MMU 



Number of EBRs 


59 


59 


109 


111 


74 


77 


Minimum length (Kbp) 


11.6 


16.9 


0.003 


0.007 


0.003 


0,007 


Maximum length (Kbp) 


4,614.2 


6,114.9 


4,505.9 


5125 


508.7 


420.3 


Median length (Kbp) 


259.3 


163.6 


51.3 


26 


51.6 


26 
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Figure 1 Experimental design of thie study. 



Recent studies in great apes have found differences in 
the genome-wide distribution of TR among species, sug- 
gesting that they might be correlated with the evolution- 
ary history of each primate chromosome [26]. More 
specifically, qualitative comparisons of TR distribution 
in great apes indicated that the TR landscape might have 
been conserved in coUinear chromosomes, but altered in 
those reorganized chromosomes [26]. Under this as- 
sumption, genomic regions that have suffered more re- 
arrangements during their evolution are expected to 
concentrate more repetitive sequences than are conserved 
regions. In fact, our dataset supports this assumption since 
rearranged chromosomes in the macaque lineage (MMU3 
and MMU5) have more TRs than do those that main- 
tained the ancestral form (MMU2, MMU7, MMUlO, 
MMU12, MMU13 and MMU18) (Mann- Whitney U 
test, p-value < 0.0001). Despite the limitations of the 
current rhesus macaque draft genome assembly and an- 
notation [47,48], this view is consistent with the lack of 
differences found in TR density between EBRs and HSBs 
in the macaque, which has maintained an ancestral karyo- 
type within Catarrhini and, consequently, can be consid- 
ered to have retained a more conserved chromosome 
complement (i.e., low degree of genome reshuffling) than 
have those of great apes [25,26]. 

Defense-responsive genes are over-represented in EBRs 

Once the evolutionary genomic landscape of the macaque 
was established, the genome-wide distribution of genes 
was further examined, paying special attention to gene 
ontology. A total of 28,595 genes was included in the ana- 
lysis: 21,023 protein-coding genes; 5,913 non-coding RNA 
genes and 1,659 pseudogenes. We scrutinized each 



macaque chromosome's complete sequence using non- 
overlapping windows of 100 Kbp in order to analyze the 
distribution of genes genome-wide. The mean distribution 
is 0.98 genes per 100 Kbp (including protein-coding genes, 
non-coding RNA genes and pseudogenes) in the whole 
macaque genome. When analyzing the distribution of 
protein-coding genes in more detail, a higher gene density 
in EBRs was detected (1.48 genes/100 Kbp), when com- 
pared with HSBs (0.73 genes/100 Kbp) (Kruskal-Wallis 
test, p-value < 0.0001). Therefore, our results indicate the 
presence of EBRs in gene-rich regions, in line with previ- 
ous observations in mammals using multi-species com- 
parative maps [22,32]. In trying to understand the reasons 
behind this pattern, initial studies reported the intergenic 
location of mammalian EBRs [22], while recent studies 
have paid special attention to the adaptive role of EBRs 
[33]. This has been the case of the pig, for example, where 
EBRs have been found to be especially rich in taste per- 
ception networks [49], suggesting that genome reshuffling 
significantly contributed to adaptation and the develop- 
ment of lineage-specific traits. Moreover, it has also been 
suggested that inversions can suppress recombination 
within the affected zones [5-8]. Consequently, rearranged 
sequences could accumulate alleles, which might be adap- 
tive for the population, and this could generate reproduct- 
ive isolation leading, eventually, to speciation. 

With this in mind, the function of the genes located at or 
near EBRs was analyzed in the macaque genome through 
the analysis of GO terms enrichment. We focused our 
study on the genes located in the positions of the EBRs de- 
tected in silico, considering flanking regions of ±200 Kbp 
in size in order to cover genes that overlapped EBR starts 
and/or ends. We applied the Functional Annotation 
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Clustering Tool from DAVID database [50] using a back- 
ground list containing 19,794 protein-coding genes. This 
was created by filtering the total protein-coding genes list 
available in the database, discarding 2,111 genes located in 
the centromeres and telomeres. Using this approach we 
found that macaque EBRs are enriched in genes related to 
the immune system. In fact, we detected a single functional 
module significandy associated to EBRs (Enrichment 
Score = 2.81, Table 2 and Additional file 4: Table S3). This 
functional module included 17 genes implicated in the 
immune response (including several GO terms, such as 
chemokine, defensin precursors and Toll-like receptor 
signalling pathway, among others) located in seven differ- 
ent macaque EBRs (Table 2, Additional file 4: Table S3). 
Interestingly, five of these genes are beta-defensins, and are 
clustered in MMUlO EBR involved in the pericentric inver- 
sion that occurred after macaque and human divergence 
(Table 2). Beta-defensins are antimicrobial peptides in- 
volved in the resistance to microbial colonization of the 
epithelial surface [51]. This cluster, whose expression is 
restricted to the male reproductive tract, is the result of a 
series of duplication events subsequently shaped by the ac- 
tion of positive selection [52,53]. In fact, previous studies 
have shown a marked differential expression of DEFB118 
and DEFB122 (two of the genes detected in MMUlO, 
Table 2) in human and macaque reproductive tissues [50]. 
Since chromosomal rearrangements are sources of genome 
variation, it might be possible that they could have influ- 
enced the structure of regulatory regions in those genes 



near EBRs. Moreover, and in line with these observations, 
we also identified an alfa-defensin cluster in the vicinities 
of the EBR detected in MMU8 (Table 2). This gene cluster 
has been previously described as an example of rapid 
evolution in primates [54]. Developing protection to new 
microbial infections is among one of the great challenges 
that species confront during the adaptation to new 
ecological niches. Thus, the generation of new defensin 
variants through different mechanisms, such as sequence 
modifications when a chromosomal reorganization occurs, 
could lead to adaptation to new environments. In fact, 
defensin clusters have also been described as being located 
in EBRs in Cetartiodactyla [55], suggesting the importance 
of genome reshuffling as an important source of new gene 
variants. 

Moreover, we also investigated genes located within 
the inverted regions that are specific for the macaque 
linage, as they could have played an adaptive role in this 
species. Statistically significant enrichment was detected 
(Enrichment Score > 1.5) for defense-response genes 
within the paracentric inversion affecting MMU5, a 
chromosomal reorganization that has been fixed in the 
lineage leading to macaque (Table 3 and Additional file 5: 
Table S4). MMU5 inversion spans approximately 33 Mbp 
across the centromere, involving one breakpoint located 
between 44.38-44.44 Mbp and the other one located be- 
tween 77.5-77.9 Mbp (Table 3). Specifically, we identified 
two different statistically significant functional modules 
consisting of 13 and 6 genes, respectively - mostiy coding 



Table 2 Genes located in macaque EBRs 


Gene 


Distance from EBR edges (Kbp) 


EBR position (chr: start-end) 


Function 


GNRHR2 


147.5 


MMUl: 126,278,886-126,466,180* 


Gonadotropin receptor 


ADIPOQ 


156.7 


MMU2: 179,026,672-179,049,651* 


Adiponectin precursor 


F7D492 


133.6 


MMU8: 8,165,750-8,338,584 


Defensin precursors 


DEFA4 


164.6 






LOC574310 


172.8 






MSMB 


31.1 


MMU9: 46,613,935-46,645,068 


Beta-microsemino protein 


DEFB123 


189.6 


MMUlO: 33,312,300-33,313,023* 


Beta-defensins 


DEFB118 


1225 






DEFB119 


113.7 






DEFB121 


152.7 






DEFB122 


1715 






CCL18 


129.4 


MMU16: 31,482,483-31,663,199* 


C-C motif chemokine precursors 


CCL23 


193.7 






CCL3 


1064 






Q76LL8 


182.8 






CCL25 


106.1 


MMUl 9: 7,76,0383-7,770,483 


C-C motif chemokine precursor 


DC-SIGN 


53 




Pathogen-recognition receptor 



Protein-coding genes of the single functional annotation ternn (Enrichment Score = 2.81 ) detected in the proximity of EBRs in the macaque genome. 
*macro-reorganizations; MMU - macaque chromosome. 
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Table 3 Genes located in the IVIMU5 inversion 



Gene Distance from EBR EBR position (start-end) 



Term 1 (ER = 


4.72) 








CXCL13 




7,5734 


44,386,92f 


i-44,442,153 


CXCL1 1 




9,135.6 






CXCLIO 




9,176.8 






CXCL9 




9,193.0 






CXCL3 




11,230.2 






PPBP 




11,276.2 






PF4 




11,282.6 






CXCLl 




11,420 






PF4V1 




11,439.3 






P67813 




11,555.7 






ALB 




11,867.3 






ODAM 




15,090.7 






KIT 




2,815.7 


77,560,1 9f 


i-77,980,505 


Term 2 (ER = 


3.52) 








SULT1E1 




15,330.9 


44,386,92f 


i-44,442,153 


UGT2B4 




17,318.6 


77,560,1 9f 


i-77,980,505 


UGT2A3 




1 7,222.2 






UGT2B33 




1 7,078.8 






UGT2B15 




1 6,769.3 






SRD5A3 




3,533.6 







Protein-coding genes contained in the two GO terms detected in thie 
macaque-specific pericentric inversion affecting chromosome MMU5. 
ER = Enrichment Score. 



for chemokines and the UDP-Glucuronosyltransferases 2B 
gene family -out of a total of 181 protein-coding genes 
present within the rearranged region (Enrichment scores: 
4.72 and 3.52, respectively) (Table 3 and Additional file 5: 
Table S4). Both functional modules detected include 
immune-response genes that might have been influenced 
by the chromosomal rearrangement in the rhesus monkey. 
Chemokines are known to play a role in the neuroinflam- 
mation process in response to infection, are present in the 
central nervous system (CNS) and are expressed in 
neurons and glial cells [56]. Moreover, there is evidence 
suggesting that they are also involved in neurodevelopment 
and neurophysiological signaling [56,57]. UDP-Glucuro 
nosyltransferases, on the other hand, are enzymes from the 
major pathway for the elimination of xenobiotics and en- 
dobiotics, and it has been suggested to play a role during 
speciation in the lineage leading to macaque [58]. 

Our observation of an over-representation of defense- 
responsive genes in both EBRs and macaque-specific in- 
versions might suggest an adaptive role of reorganizations 
in this species. Previous studies have reported that a small 
proportion of the mammalian genome, i.e., 4% in the case 
of the human genome, is under selective constraints, espe- 
cially so for coding regions, introns and intergenic regions 



[59]. This suggests that in certain regions the fitness cost 
is so pronounced (i.e., could be lethal or deleterious for 
the individual and the progeny) that rearrangements are 
not allowed (i.e., [60]). But it has also been shown that this 
constraint could be somewhat relaxed in the promoters of 
genes linked to the immune system, reproduction and 
perception [59], allowing for the generation of new vari- 
ability to ensure adaptation to new environments. In light 
of our results, this might be the case for the genomic 
regions under study in the macaque. However, whether 
immune-related genes are directly involved in lineage- 
specific adaptation, as has been previously suggested for 
macaque [58], needs further validation. 

Genome reshuffling and its effect on chromosome- 
specific recombination landscapes 

Once the genomic structural changes were defined in 
the rhesus macaque, together with the genome distribu- 
tion of coding-genes across evolutionary regions, we fur- 
ther experimentally analyzed the meiotic recombination 
landscape as a way to explore the adaptive role of chromo- 
somal changes. Under the "suppressed recombination" 
model of chromosomal evolution, chromosome rearrange- 
ments would have a minimal influence on fitness, but 
would rather suppress recombination within the genomic 
regions affected, thus contributing to the accumulation of 
gene incompatibilities [2,5-9]. The analysis of the most re- 
cent human and chimpanzee recombination maps inferred 
from genome-wide single-nucleotide polymorphism (SNP) 
data revealed that the standardized recombination rate was 
significantly lower in rearranged rather than in coUinear 
chromosomes [5]. In the case of rhesus macaque, chro 
mosome-specific recombination maps are available for very 
few chromosomes [61], and whether or not chromosomal 
reorganizations that have been fixed in the macaque 
lineage have affected the recombination landscape was ad- 
dressed in our study. 

Catarrhini monkeys, and more specifically the Tribe 
Papionini, are characterized by the presence of highly con- 
served karyotypes in terms of diploid number and chromo- 
some homologies [62,63]. It has been described that 
species sharing the same karyotype (i.e., no major genome 
reshuffling) maintain the chromosomal distribution of mei- 
otic crossover in conserved chromosomes on a broad scale 
(Mbp, the resolution provided by the in situ immunoflour- 
escence detection of meiotic proteins) [64,65]. Therefore, 
we expect that the recombination pattern in Macaca and 
Cercocebus chromosomes (species that share the same 
karyotype and belong to the Tribe Papionini [62,63]) are 
conserved at the Mbp resolution, thus allowing us to ex- 
trapolate the results obtained in both species, as previously 
described [64,65]. In this context, and taking advantage of 
the EBRs detected in silico in the macaque genome, we 
subsequendy studied the chromosomal distribution of 
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meiotic crossover (COs - here exemplified as MLHl foci) 
in macaque chromosomes affected by inversions since their 
Catarrihini common ancestor (MMU5). This was done by 
combining the immunostaining of meiotic proteins and 
fluorescence in situ hybridization (FISH) with specific BAC 
clones spanning the EBRs detected in our in silico scanning 
(Figure 2, Table 4 and Additional file 6: Table S5). The ana- 
lysis of COs distribution according to centromere positions 
and EBR-specific BAC probes allowed us to experimentally 
determine the chromosome position of evolutionarUy reor- 
ganized regions directly on pachytene chromosomes along 
chromosomal axes. Our efforts were concentrated on three 
chromosomes, which represented three different evolution- 
ary states: (i) a chromosome with a pericentric inversion 
(CT05/MMU5) specific for the Cercocebus and Macaca 
lineages; (ii) a paracentric inversion (CT09/MMU9) that 



has been fixed in human and chimpanzee lineages, and 
therefore has been maintained coUinear in the macaque 
lineage; and (iii) a coUinear chromosome (CT06/MMU6) 
that has been maintained coUinear since the Catarrhini 
common ancestor. With this in mind, our main goal was 
two-fold: (i) reconstruct chromosome-specific recombin- 
ation maps, and (ii) encompass the distribution of meiotic 
crossovers across the inverted regions in order to test 
whether chromosomal inversions have an effect in reducing 
recombination. 

A total of 258 spermatocytes at pachynema was ana- 
lyzed in order to obtain chromosome-specific recombin- 
ation maps. We detected an average of 35.24 (±3.01) 
MLHl foci per analyzed cell (Figure 2). Assuming that one 
MLHl focus can be translated into 50 centimorgans (cM) 
in genetic length according to previous studies [61,64-66], 
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Figure 2 Recombination features in Cercocebus torquatus. (A) Number of MLHl foci detected per cell in Cercocebus torquatus (CTO). The 
horizontal bar indicates the mean. (B) Comparison of the percentage of cells with a different number of MLHl foci (crossovers, CO) in the 
different chromosomes analyzed: CT05, CT06 and CT09. (C-D) Sequential image of a spermatocyte at pachynema from Cercocebus torquatus 
depicting (C) triple immunostaining against SYCP3 (red), MLHl (green) and centromeres (blue). (D) Representation of the same cell after applying 
sequential fluorescent in situ hybridization (FISH) with the BAC probe RPl 1-926D9, specific for CT09. 
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Table 4 Chromosome-specific recombination analysis 



Chr 


N 


(jm SC 


MLH1 
foci/ 
SC 


MLH1 

foci/ 
pm 


MLH1 foc//pm 


p-arm 


q-arm 


Inside inv 


Outside inv 


CT05 


87 


3.85 


4.01 


2.01 


0.26 


0.04 


031 


CT05p 


85 


2.25 


6.12 


2.25 


0.27 


0.06 


030 


CT06q 


85 


2.26 


6.12 


2.25 


0.27 


0.09* 


036 


CT09 


86 


3.19 


5.56 


1.99 


0.23 


0.03* 


031 



Chr: Chromosome, CTO: Cercocebus torquatus, N: number of cells analyzed, SC: 
synaptonemal complex, Inv: inverted region. CT06p: inversion simulated In the 
small chromosome arm. CT06q: inversion simulated In the long chromosome 
arm. # CT06 density {expressed as MLH1 foc//pm) Inside and outside the 
inversion refers to the simulated region for the comparison in each case. 
*indicates statistical significance (p-value < 0.05) when comparing reorganized 
region with simulated inversion in the colllnear chromosome CT06 (Mann 
Whitney U test). 

these recombination events correspond to an autosomal 
genetic length of 1,760 cM for C. torquatus, very similar 
to the 1,950 cM obtained in M. mulatta [61]. Our results 
also indicated that all chromosome pairs present, at least, 
one recombination event, with an average number of 1.8 
COs per homologous autosomic chromosome and 0.88 
COs per autosomic arm. Moreover, a maximum of three 
recombination events in the large chromosomes was 
detected (Figure 2). In addition, for each chromosome an- 
alyzed, its total length (expressed in |im) and the density 
of COs (MLHl focil \im) per chromosome and chromo- 
somal arm were calculated (Table 4). All three chromo- 
somes presented equivalent CO densities along their 
synaptonemal complexes: 0.26 COs/|im for CT05, 0.27 
COs/^m for CT06 and 0.23 COs/^tm for CT09. In 
general, we detected an increase in CO frequency to- 
wards the chromosomal telomeric regions, whereas 
centromere regions showed very low CO frequencies 
(Figure 3), mirroring previous studies in mammalian 
species [61,64-66]. 

Subsequendy, CO density was analyzed inside and out- 
side each inverted region using the selected BAC probes to 
label the location of the breakpoint in each chromosomal 
region affected by the inversion (Figures 2 and 3). As a 
general trend, we observed low CO densities within 
inverted regions (0.04 and 0.03 MLHl focil\i.m for CT05 
and CT09, respectively), when compared to chromosomal 
regions outside the reorganized area (0.31 COs MLHl foci 
per |im in both cases) in both rearranged chromosomes 
(Table 4). When analyzing the recombination rate within 
inverted regions among chromosomes, significant reduc- 
tion of CO density was observed in CT05 and CT09 
within the inverted regions, when compared with regions 
outside the inversion (Table 4, Mann- Whitney U test, p- 
value < 0.05). Such differences were not observed in the 
simulated inversion in CT06, a chromosome that has been 
maintained coUinear in the macaque lineage (Table 4). 

Subsequendy, we tested whether the suppression of re- 
combination observed within reorganized areas was due to 



the low recombination rate characteristic of pericentromeric 
regions. To do so, we considered as a pericentromeric 
region an area extending 30% of each chromosome arm 
from the centromere towards the telomeric region and 
compared it to the CO density observed in the same re- 
gion in the collinear chromosome (CT06). When all 
chromosomes were compared, our results showed no stat- 
istical differences among the pericentromeric regions in 
the small arms (CT05p, CT06p and CT09p) (Kruskal- 
Wallis test, p-value > 0.05); however, differences in CO 
density were significant when considering the long arms 
(CT05q, CT06q and CT09q) (Kruskal-Wallis test, 
p-value = 0.018). Moreover, a significant reduction was ob- 
served in the CO density within the rearranged region in 
CT09, when compared to the collinear chromosome 
CT06 (Mann- Whitney U test, p-value = 0.016), (Figure 3 
and Table 4). These differences were not significant, how- 
ever, when the inverted region of CT05 was compared to 
the collinear chromosome CT06 (Mann- Whitney U test, 
p-value > 0.05). 

What are the evolutionary implications of our observa- 
tions in light of the "recombination suppression" model? 
Despite the fact that the "recombination suppression" 
was initially proposed to explain differences in recom- 
bination rates within reorganized genomics in heterokar- 
yotypes (i.e., heterozygotes), we observed a reduction of 
recombination in fixed rearrangements, raising intri- 
guing questions about the mechanisms involved. Previ- 
ous studies in great apes have revealed that rearranged 
chromosomes presented significantly lower recombin- 
ation rates than do chromosomes that have been main- 
tained collinear since a common ancestor, and this was 
related to the lineage in which they become fixed [5]. 
Importantly, inverted regions had lower recombination 
rates than did collinear and non-inverted regions, inde- 
pendently of the effect of centromeres [5]. Although at 
this stage it would be premature to argue about the evo- 
lutionary forces behind this pattern, our results highlight 
the importance of the study of recombination framed by 
the evolutionary history of chromosomes and, in greater 
extent, genomes. Incorporating more chromosomes into 
the experimental study would be necessary to detect a 
clearer genomic effect of the inversions in the distribu- 
tion of recombination patterns. 

Conclusions 

Genomic rearrangements might play an important role in 
local adaptation and species divergence by the modifica- 
tion of both the structure and regulation of genes located 
near the affected regions. Here, we provide a highly 
refined description of the chromosomal reorganizations 
and evolutionary breakpoint regions in the human and 
rhesus macaque genomes based on orthologous genes and 
genome sequence alignments. The high-resolution map of 
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CTO (left) and HSA (right): RPl 1-779 N22 and RP1 1-8 N8 for CT05, RPll-157 FT for CT06 and RPl 1-926D9 for CT09. 



EBRs and HSBs defined in the macaque genome has 
revealed an interesting pattern: evolutionary break- 
points are gene-rich regions, with a significant func- 
tional clustering for genes related to the immune 
system. Furthermore, and in light of our observations, 
inversions can induce a reduction in the recombination 
rate among the different alleles contained in the inver- 
sion, which could be favoring adaptive alleles. Future 
comparative research on the effect of chromosomal 
reorganization on recombination as outlined above 
should be an effective means to enhance our knowledge 
of the role of genome reshuffling in evolution. 



Methods 

Whole-genome comparisons and evolutionary 
breakpoint definition 

The experimental design of the study is represented in 
Figure 1. In order to detect the evolutionary breakpoint 
regions (EBRs) and homologous synteny blocks (HSBs) 
between human and macaque whole-genome sequences, 
two different algorithms were applied: SyntenyTracker 
[67] and Cassis [68]. Both approaches compare the order 
and orientation of orthologous markers (genes) between 
genomes, detecting changes both in the sequence order 
and the position of the HSBs and EBRs. Orthologous 
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genes between human (GRch37.p7) and rhesus monkey 
(MMUL_1.0) genomes were obtained through the 
BioMart database (http://www.ensembl.org/index.html). 
SyntenyTmcker [67] determines the position in the 
chromosome sequences of both genomes, providing infor- 
mation of the relative orientation of each HSB. This en- 
ables the detection of chromosomal rearrangements such 
as inversions and/or translocations between two genomes. 
Once the HSBs were detected, genomic regions between 
consecutive HSBs were considered EBRs. Cassis [68], on 
the other hand, is especially designed to define breakpoint 
regions, providing information about different types of re- 
arrangements, such as inversions, translocations or indels. 
Both algorithms were applied as previously described [5] 
using default parameters for SyntenyTmcker and level 1 in 
the lastz alignment in Cassis. In order to obtain the gen- 
omic positions of EBRs in both genomes, the analysis was 
performed in two directions: (i) using both the human 
genome as reference genome, and (ii) the macaque 
genome as reference. Following previous studies [19,22], 
we considered EBRs that were 4 Mbp in size or less. Re- 
gions larger than 4 Mbp in size were considered "gaps". 
Furthermore, we labeled as telomeric/subtelomeric the 2 
Mbp at the ends of each human chromosome and as 
centromeric/pericentromeric the 2 Mbp regions flanking 
the unknown nucleotides (Ns), as previously described [25]. 

Once the genomic positions of EBRs were obtained, we 
followed conservative criteria in order to avoid false posi- 
tives. To do so, EBRs located at telomeres, centromeres 
and gaps were excluded from the analysis. The resulting 
EBRs were classified according to whether they are in- 
volved in macro-rearrangements (rearranged regions > 1.4 
Mbp) or micro-rearrangements (rearranged regions < 1.4 
Mbp). Simultaneously, we also classified each EBR de- 
pending on which type of chromosomal rearrangement 
was involved, that is, inversion, fusion or fission, following 
previous studies [30]. 

Analysis of repetitive elements and gene screening 

The distribution of TR in the macaque genome using the 
eTandem algorithm was analyzed (part of EMBOSS 6.0.1 
[69]). The eTandem algorithm was run with a minimum 
repeat unit of 2 bp and a maximum repeat unit of 100 bp, 
as previously described [26]. The resulting output files 
were computed for the detection of overlapping TR, and 
the canonical motif was defined using in-home Perl 
scripts. In order to compare the distribution of TR along 
macaque chromosomes, we counted the base pairs of tan- 
dem repeats in 100 Kbp windows for each chromosome. 
Finally, each window was labeled according to its position: 
telomere, centromere, HSBs or EBRs. Using Perl scripts, 
we computed the density of TR and merged the positions 
of TR with the different types of genomic regions in the 
human genome. 



Then, the number and the genomic position of anno- 
tated genes in the macaque genome were considered 
(RefSeq from the MMUL l.O assembly) to perform the 
gene distribution analysis using the BioMart browser of 
Ensembl (release 67). We grouped all genes with a known 
function in non-overlapping windows of 100 Kbp and 
labeled each window according to its position: telomere, 
centromere, HSBs or EBRs. In this case, the variable 
considered was gene count in each window, given that we 
analyzed presence/absence of genes, rather than the num- 
ber of base-pairs covered by genes. 

The Functional Annotation Clustering tool in DAVID 
{Database for Annotation, Visualization, and Integrated 
Discovery, v6.7) [50], was used in order to identify over- 
represented biological terms contained in EBRs. Func- 
tional annotation clustering allows for the biological 
interpretation at a 'biological module' level of the most 
relevant biological terms (GO). Following algorithm's rec- 
ommendations, all clusters analyzed included a minimum 
of 10 genes and a maximum of 3,000 [50]. In DAVID an- 
notation system, Fisher Exact is adopted to measure the 
gene-enrichment in annotation terms by means of an 
EASE-score, a modified Fisher Exact P-Value. EASE- 
scores equal or smaller than 0.05 were considered statisti- 
cally significant (i.e., strongly enriched in the annotation 
categories). Additionally, the system uses the group En- 
richment Score (a geometric mean of member's p-values 
in a corresponding annotation cluster) to rank the bio- 
logical significance of the genes found in a cluster. Enrich- 
ment Scores > 1.5 indicated significant over-represented of 
gene functions. 

Biological samples 

Metaphase chromosomes were prepared from peripheral 
blood samples obtained from one female rhesus macaque 
(Pare Zoologic de Barcelona, Spain). Cultures from per- 
ipheral blood samples were processed under standard 
conditions in order to obtain chromosome preparations as 
previously described [35]. Additionally, testicular tissue 
from an adult individual of Cercocebus torquatus (CTO, 
2n = 42) with proven fertility was used for the study of 
meiotic recombination. In order to obtain spermatocyte 
spreads, testicular tissue was processed as previously de- 
scribed [64,65]. 

Immunofluorescence 

Immunostaining of meiocytes was performed as previously 
described [64,65]. Different sets of antibodies were used: 
rabbit anti-SYCP3 (Abeam), human anti-CenP (human 
serum CREST, a kind gift from Dr. M. Fritzel) and mouse 
anti-MLHl (Pharmigen) for the detection of meiotic cross- 
overs (COs), all of them diluted in PTBG solution (0.05% 
Tween 20 in PBS) 1:200, 1:200 and 1:100, respectively. 
Fluorochrome-conjugated secondary antibodies (all from 
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Jackson Immunoresearch) were used for detection: goat 
anti-rabbit conjugated with Cy3, goat anti-human conju- 
gated with Cy5 and goat anti-mouse conjugated with FITC 
diluted 1:100 in PTBG. 

Fluorescence in situ hybridization 

BAG clones spanning evolutionary breakpoints were ob- 
tained from the human library available at CHORI (Chil- 
dren's Hospital Oakland Research Institute) (Additional 
file 6: Table S5). DNA from BACs was extracted according 
to standard protocols using a commercial kit (QIAGEN 
Plasmid). Fluorescence in situ hybridization (FISH) with 
specific BAG clones was performed on both metaphase 
chromosomes and spermatocyte spreads as previously 
described [40,64]. Briefly, 1 |ig of the DNA plasmid was 
labeled with dUTP-digoxygenine by Nick Translation 
(Abbot kit) and ethanol precipitated with competitor 
DNA (Got-1 human DNA, Invitrogen, 1 mg/ml), salmon 
sperm DNA (Invitrogen, 10 mg/ml) and 1/10 volume of 
3 mol/L sodium acetate overnight at -20°G. The precipi- 
tated probe mix was resuspended in 14 ml hybridization 
buffer (50% deionized formamide, 10% dextran sulfate, 
2xSSG and 0.5 mol/L phosphate), denatured 80°G for 
10 min and pre-annealed at 37°G for 1 h. Preparations 
were visualized using a Zeiss Axioskop epifluorescence 
microscope equipped with the appropriate filters and 
a charged coupled-device camera (ProgRes" GSlOplus, Je 
noptik). 

Recombination analysis 

The Micromeasure 3.3 software [70] was used for image 
analysis and for the construction of chromosome-specific 
recombination maps based on the relative distances be- 
tween adjacent MLHl foci, a marker for meiotic cross- 
overs (GOs). For each chromosome analyzed, the position 
of each MLHl foci was recorded as a relative position (as 
the percentage of total length of the synaptonemal com- 
plex, SG) from the centromere, identified by the centro- 
meric signal in each preparation as described previously 
[64,65]. Using the centromere as a reference, the positions 
of each MLHl focus were calculated along the SG, from 
the centromere to the telomere. Thus, for comparison 
among chromosomes, the position of MLHl foci was 
expressed as the relative position of each GO to the length 
of the chromosome (the length of each SG was divided 
into 10% intervals). To convert the MLHl foci to genetic 
distances, the number of MLHl foci detected per SG was 
multiplied by a factor of 50 map units (1 crossover = 
50 cM) [64,65]. 

The effect of chromosome inversions on the GO distri- 
bution pattern was analyzed by calculating MLHl foci 
density within inverted and non-inverted regions, consid- 
ering the length (expressed in \im) for each region, so the 
differences due to SG lengths for each chromosome were 



normalized. In order to delimit the inversions directly in 
spermatocytes, the centromere position and the specific 
BAG probes labeling the breakpoint distal to centromere 
was used (Figure 1 and Additional file 6: Table S5). To 
allow for comparison among chromosomes, GO distribu- 
tion was expressed as the relative position of each GO to 
the length of the chromosome (the length of each SG was 
divided into 10% intervals). 

Furthermore, and in order to disentangle the centro- 
meric effect on recombination, the GO distribution was 
compared between the rearranged chromosomes and the 
collinear one, using the last as a control. To this purpose, 
we simulated an inversion in the collinear chromosome, 
by analyzing the recombination pattern of a region of 
(proportionally) the same size as the observed inverted re- 
gion. We constructed plots of cumulative frequency to 
study the pattern along the chromosome arms. 

Statistical analysis was performed using JMP 10 soft- 
ware (SAS Institute Inc.) and IBM SPSS Statistics 20, 
using the Kolmogorov-Smirnov-Lilliefors test for normal- 
ity, and the Kruskal-Wallis and Mann-Whitney U tests 
for comparisons. 

Additional files 



Abbreviations 

HSBs: Homologous syntenic blocks; EBRs: Evolutionary breakpoint regions; 
SYCP3: Synaptonemal complex protein 3; SC: Synaptonemal complex; 
MLHl: MutL homolog 1; BAC: Bacterial artificial clones; CHORI: Children's 
Hospital Oakland Research Institute; COs: Crossovers; MMU: Macaca mulatto; 
CIO: Cercocebus torquatus; TR: Tandem repeats; HSA: Homo sapiens. 



Additional file 1: Table SI. EBR positions involved in the macro- 
rearrangements (inversions spanning more than 4 Mbp, fusions and 
fissions) between human and macaque genomes detected in our study. 

Additional file 2: Table S2. EBR positions involved in the micro- 
rearrangements (inversions spanning less than 4 Mbp, indels and high- 
complex regions) between human and macaque genomes detected in 
our study. 

Additional file 3: Figure SI. Highly refined map of the reorganizations 
and evolutionary breakpoint regions in the human and rhesus macaque 
genomes. Representation of HSB (gray blocks) and EBRS (white regions) 
between human and rhesus monkey, using human as the reference 
genome detected by SyntenyTracker (Macaque_Synteny) and Cassis 
(IVlacaque_Cassis) algorithms, as well as the final model (Macaque_Einal) 
(from left to right in each chromosome representation). The final model 
is the result of merging the outputs of both programs. Inset numbers 
represent the homologous rhesus monkey chromosomes. Hatched areas 
represent heterochromatin in the human genome. 

Additional file 4: Table S3. Output from Functional Annotation Tool of 
DAVID using genes inside EBR positions (H — 200 kbp) (see materials and 
methods for details). Annotation clusters are statistically significant when 
Enrichment Score is above 1.5. 

Additional file 5: Table S4. Output from Functional Annotation Tool of 
DAVID using genes inside the IV1MU5 (see materials and methods for 
details). Annotation clusters are statistically significant when Enrichment 
Score is above 1 .5. 

Additional file 6: Table S5. Selection of human BAC clones for the 
EBRs detected in silico. 
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